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ABSTRACT 


An  Interactive  computer  model  of  a  highly  enriched  pressurized  water 
reactor  was  developed,  using  the  applicable  plant  parameters  from  the 
Shlpplngport  Atomic  Power  Station.  The  point  reactor  kinetics  equations 
for  one  delayed  neutron  precursor  group  were  linearized  using  small 
perturbation  theory.  The  model  Included  both  moderator  and  Xenon-135 
reactivity  feedback  effects,  as  well  as  an  automatic  reactor  protection 
and  average  reactor  coolant  temperature  control  system.  The  thermal 
response  of  the  model  plant  was  simulated  for  normal  operating  transients 
Induced  either  by  control  rod  or  turbine  load  changes.  The  post  shutdown 
Xenon  transient  response  was  also  modeled.  The  interactive  program  was 
coded  In  FORTRAN- IV  language,  and  the  simulation  program  was  coded  In 
IBM  CSMP-III  language. 
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I.  INTRODUCTION 


Computers  have  been  used  In  the  design  and  analysis  of  nuclear 
reactors  since  the  Inception  of  reactor  technology,  and  many  codes  have 
been  developed.  More  recently,  digital  computer  programs  have  been  used 
as  learning  devices  for  nuclear  engineering  students.  The  formulation 
of  the  computational  problem  of  predicting  the  behavior  of  a  nuclear  re¬ 
actor  has  been  greatly  facilitated  with  the  use  of  high  level  computer 
languages. 

The  purpose  of  this  work  was  to  develop  a  computer  assisted  learning 
device  to  be  used  by  students  taking  nuclear  engineering  courses  at  the 
Naval  Postgraduate  School.  The  program  graphically  displays  the  simulated 
kinetic  and  thermal  transient  responses  of  a  pressurized  water  reactor 
power  plant.  The  reactivity  feedback  effects  of  Xenon-135  poisoning  and 
moderator  temperature  are  Incorporated  into  the  model.  Normal  operating 
transients,  starting  from  a  steady  state  critical  condition,  can  be 
induced  by  ramp  changes  in  either  control  rod  position  or  turbine  load. 

The  initial  reactor  power  level  and  reactivity  change  mechanism  are 
chosen  by  the  program  user.  For  user  convenience,  these  inputs  are 
prompted  and  entered  interactively,  after  which  the  simulation  is  run 
with  no  further  user  action  required. 

A  pressurized  water  reactor  (PWR)  was  modeled  as  it  is  the  most 
common  type  of  reactor  used  for  power  plant  application.  As  a  realistic 
reference,  the  model  design  incorporated  the  applicable  characteristics 
of  the  Shippingport  Atomic  Power  Station,  whose  nuclear  parameters  were 
the  most  compatible,  of  available  PWR  core  data,  with  the  assumption  of 
a  highly  enriched  core  made  in  the  model  design. 
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The  model  was  developed  by  treating  the  reactor  kinetics  and  other 
plant  component  thermodynamic  relationships  as  transfer  functions  [1]. 
The  point  reactor  kinetics  equations  were  linearized  using  a  Taylor  ex¬ 
pansion  or  small  perturbation  technique  [2].  This  same  perturbation 
technique  and  lumped  parameter  analysis  were  applied  to  the  plant's 
linear  differential  heat  transfer  and  Xenon-135  equations  [3].  The 
individual  transfer  functions  were  then  Integrated  into  an  overall  plant 
block  diagram.  The  program  also  features  a  reactor  protection  and 
average  reactor  coolant  temperature  control  scheme. 

IBM  CSMP-III  language  was  chosen  to  formulate  the  model  simulation 
as  It  has  Inherent  routines  to  Invert  the  transfer  functions  to  the  time 
doma  i  n . 


I 
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II.  MODEL  DESIGN  CONSIDERATIONS 


Lumped  parameter  analysis  and  first  order  perturbation  theory  were 
used  throughout  the  model  development.  All  the  nuclear  variables  which 
were  used  in  the  reactor  kinetics  and  Xenon-135  decay  equations  were 
considered  to  be  averaged  values  of  the  variables  over  the  neutron 
energy  spectrum.  Similarly,  the  thermodynamic  variables  which  were 
used  in  the  plant's  heat  transfer  equations  were  considered  to  be 
averaged  values  over  the  volume  of  the  particular  component. 

Because  of  these  approximations  and  other  assumptions  made  in  the 
model  development,  the  simulation  should  not  be  considered  a  design  or 
stability  analysis.  Instead,  the  simulation  shows  the  model's  large 
scale  reactor  kinetic  and  thermal  trends  during  normal  operating 
transients. 

For  similar  reasons,  although  certain  plant  parameters  of  the 
Shippingport  Atomic  Power  Station  were  incorporated  into  the  model,  the 
simulation  cannot  be  considered  to  reflect  the  operating  characteristics 
of  this  power  plant.  The, Shippingport  core  is  a  seed  and  blanket  type 
whereas  the  model  core  is  a  uniform  mixture  of  highly  enriched  fuel 
and  support  materials. 

A.  PRESSURIZED  WATER  REACTOR 

The  pressurized  water  reactor  (PWR)  is  the  most  widely  used  reactor 
type  in  central  power  plant  applications  and  the  only  type  currently 
used  in  naval  propulsion. 

A  simplified  schematic  of  the  modeled  PWR  plant  is  shown  in  Figure 
1.  The  modeled  reactor  contains  a  highly  enriched  core  which  is  light 
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Figure  1.  Schematic  of  the  Model  Plant 

water  moderated  and  cooled.  The  modeled  plant  contains  two  closed  loop 
thermodynamic  systems  coupled  by  a  heat  exchanger. 

In  the  primary  system,  heat  generated  from  thermal  fission  is  trans¬ 
ferred  to  the  coolant  as  it  passes  through  the  reactor,  raising  the 
coolant's  temperature.  The  high  temperature  coolant  leaving  the  reactor 
is  circulated  by  a  pump  through  tubes  inside  of  a  heat  exchanger  where 
it  gives  up  some  heat  to  a  secondary  fluid.  A  pressurizer  maintains 
the  primary  coolant  at  a  sufficiently  high  pressure  so  that  the  bulk 
temperature  of  the  coolant  is  kept  below  the  saturation  temperature. 

On  the  secondary  side  of  the  heat  exchanger,  the  temperature  of  the 
entering  feedwater  Is  raised  to  the  boiling  point  and  saturated  steam 
is  produced.  The  steam  Is  delivered  to  a  turbine,  is  condensed  and 
returned  to  the  heat  exchanger,  completing  the  second  closed  loop. 


12 


B.  MODEL  AND  PROGRAM  FLEXIBILITY 


The  model  was  designed  to  simulate  only  normal  operating  transients 
from  an  Initial  steady  state  condition.  Reactor  accidents  and  startup 
were  not  considered  In  the  model  development.  However,  the  model's  re¬ 
actor  control  module  will  simulate  either  a  full  or  constant  Insertion 
of  control  rods  if  certain  parameters  are  exceeded.  This  feature  was 
not  Incorporated  for  accident  analysis,  but  to  keep  parameters  within 
the  limits  of  the  assumptions  made  in  the  model  development. 

The  applicable  Shippingport  plant  characteristics  are  fixed  In  the 
simulation  program.  User  inputs  are  limited  to  choosing  Initial  power 
level,  from  which  other  initial  parameters  are  adjusted,  and  the  plant 
perturbation  mechanism,  either  control  rod  movement  or  turbine  load 
change.  These  perturbations  occur  at  fixed  rates  which  limit  the  amount 
of  reactivity  that  can  be  inserted  during  the  simulation. 

Originally  it  was  envisioned  that  an  overall  program  would  feature 
not  only  interactive  input  capability  but  also  produce  real  time  graphi¬ 
cal  displays  of  the  transient  for  a  user  at  a  time  sharing  computer 
terminal  which  had  graphics  capability.  However,  the  computer  language 
(IBM  CSMP-III)  required  to  solve  the  model's  algorithms  was  not  avail¬ 
able  on  the  time  sharing  system,  and  furthermore,  the  programs  long 
execution  time  precluded  any  real  time  response. 

In  order  to  still  provide  as  much  user  facility  as  possible  with 
these  restrictions,  the  interactive  feature  was  partially  retained  for 
data  input.  Using  a  computer  language  available  on  the  time  sharing 
system  (FORTRAN-IV),  a  program  which  interactively  prompts  the  user  to 


*, 
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enter  specified  controlling  parameters  was  developed.  This  program  In¬ 
corporates  logic  routines  which  permit  Inputs  only  within  requested 
ranges.  This  feature  allows  for  data  reentry  If  a  user  error  Is  made, 
and  ensures  that  only  inputs  compatible  with  the  model  development's 
assumptions  are  used  in  the  simulation. 

After  the  input  Is  completed,  the  initial  transient  conditions  are 
displayed  at  the  terminal.  A  separate  internal  control  program  then 
transfers  the  user  supplied  Inputs  into  the  simulation  program.  The 
control  program  then  transfers  the  simulation  program  to  the  batch 
processing  system  for  execution.  The  interactive  program  notifies  the 
user  that  this  has  been  done.  Hard  copy  plots  of  the  transient  are 
subsequently  produced.  The  post-shutdown  Xenon  behavior  module  is 
located  in  the  interactive  program  and  produces  real  time  graphical 
displays  at  the  terminal. 

Thus  there  are  three  distinct  programs: 

1.  An  interactive  program  to  prompt  and  receive  user  input, 
and  also  simulate  post-shutdown  Xenon  behavior. 

2.  A  batch  processed  program  which  simulates  the  transient 
response  to  the  user  inputs 

3.  An  internal  control  program  which  interfaces  these  two 
programs . 
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III.  MODEL  DEVELOPMENT 


A.  POINT  REACTOR  KINETICS  EQUATIONS 

A  lengthly  and  formal  derivation  of  the  point  kinetics  equations  is 
found  in  Reference  4  and  win  not  be  repeated  here.  More  generally  these 
equations  are  obtained  from  one  group  neutron  diffusion  theory  with  the 
assumption  that  the  neutron  flux  Is  deparable  in  time  and  space,  and 
with  the  inclusion  of  delayed  neutrons  [5].  They  are  listed  below. 

djM  .  p(t)  ,  £Nci(t)  +  Q(t)  (1) 


dCt(t) 

~dt 


AfCjU) 


i  -  1.2 


*  •  •  • 


(2) 


where  j!(t) 
/>(t) 

c1(t) 

Q(t) 

* 

and  /i 


Flux  amplitude  function 
Reactivity 

Effective  concentration  of  the  1th  delayed 
neutron  precursor  group 

Extraneous  delayed  neutron  source  strength 

Effective  delayed  neutron  fraction  from 
the  ith  group 

Decay  constant  of  the  1th  group 

ZA  3  Total  effective  delayed  neutron 
*  fraction 


These  equations  are  referred  to  as  the  point  reactor  kinetics  equa¬ 
tions  not  because  the  reactor  Is  considered  as  a  single  point  in  their 
application,  but  since  spatial  variations  are  neglected.  and 

are  assumed  to  be  constant. 

If  the  neutron  flux  spatial  variation  Is  assumed  to  be  time  Invari¬ 
ant,  0(t)  can  be  considered  to  represent  the  number  of  neutrons  in  the 


core.  Equation  (1)  may  then  be  seen  as  a  neutron  rate  equation  where 
the  three  terms  on  the  right  hand  side  represent  the  rate  of  production 
of  prompt,  delayed,  and  source  neutrons  respectively. 


Equation  (2)  Is  a  rate  equation  for  the  1th  delayed  neutron  precursor 

group.  The  two  terms  on  the  right  hand  side  represent  production  and 

decay  rates  respectively.  The  precursors  are  comprised  of  approximately 

thirty  Isotopes  which  historically  have  been  divided  Into  six  groups 

-1  235 

with  decay  constants  ranging  from  0.0124  to  3.0  seconds  for  U.  In 
the  model,  the  precursors  were  considered  to  be  represented  by  one 
effective  group  characterized  by  an  averaged  decay  constant  [2]. 


and  an  effective  delayed  neutron  fraction 

& 

In  equations  (1)  and  (2),  reactor  power  may  be  substituted  for  neutron 
flux  If  the  precursor  concentration  Is  modified  by  c  -  Ef  S  fCold 

since 

P-Ef  Zft 

where  P  ■  Reactor  power 

E^  ■  Energy  released  per  fission 

«  Macroscopic  fission  cross  section 
0  ■  Integrated  one  group  neutron  flux 

For  a  reactor  operating  at  power,  the  source  term  contribution  to 
the  overall  neutron  population  is  negligible,  l.e.,  Q  *  0. 


16 


Incorporating  these  assumptions,  equations  (1)  and  (2)  reduce  to: 

p(t)  ♦  Ac(t)  (3) 

dCitl  .  -£p(t)  -  AC(t)  (4) 


1 •  Zero  Power  Reactor  Transfer  Function 

Consider  an  Initially  critical  steady  state  reactor  at  some  power 
level  PQ  at  t  2f  0.  Equations  (3)  and  (4)  each  yield: 


Co"JOT  Po 

Now  let  P(t)  »  PQ  +  ^>(t) 
C(t)  -  CQ  +  tf fc(t) 
/»(t)  -/f0  +  <0(t) 


(5) 

(6) 


where  the  zero  subscript  denotes  the  initial  steady  state  value 
and  the  delta  prefix  a  small  perturbation  Imposed  at  t  ■  0  about 
this  value. 

Noting  that  /O  ■  0  for  a  critical  reactor  and  substituting  equa¬ 
tions  (5)  and  (6)  into  equations  (3)  and  (4)  yields  for  t  i  0 


|t  ^P(t)  -  £  *>(t)-  -f-P(t)+ 

|t  ^CCt)  -  «/Vct)  -  \</t(t) 

where  the  term  In  equation  (7)  has  been  neglected. 

Talcing  the  Laplace  transforms  of  equations  (7)  and  (8): 

s  ^>(s)  ,  JLfyW  -  JL**P (*)  4*  kJc<a\ 

% 

Sc(»)  •  ~^SP(0)  -  A</*CU» 


(7) 

(8) 


(9) 

00) 
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Solving  equation  (10)  for  C(s),  substituting  this  Into  equation 
(9),  and  rearranging,  yields: 

ft 

s(a» 

This  Is  the  zero  power  reactor  transfer  function,  so  called  because 
the  reactor  Is  assumed  to  be  operating  at  a  sufficiently  low  enough 
power  that  no  feedback  effects  are  realized.  These  effects  are  examined 
In  the  following  sections. 

Equation  (11)  may  be  represented  In  block  diagram  form  as  shown  In 
Figure  2. 


Figure  2.  Block  Diagram  of  the  Zero  Power  Reactor 
Transfer  Function 


where  Z(s) 


S(A 


2.  Limitations  of  the  Point  Reactor  Kinetics  Equations 

The  point  reactor  kinetics  equations  derivation  Is  based  on 
the  assumption  that  the  spatial  dependence  of  the  neutron  flux  Is 
negligible.  This  assumption  limits  the  validity  of  these  equations  to 
transients  where  this  remains  a  reasonable  approximation  such  as  those 
which  result  In  only  small  changes  In  reactivity.  The  small  perturba¬ 
tion  technique  used  In  the  development  of  the  zero  power  reactor 
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transfer  function  also  requires  only  small  changes  In  reactivity  If  the 
first  order  approximation  Is  to  hold.  Specifically,  must  be  less 
than  O.S/3  .  Physically,  when  /01$  greater  than  /&  the  reactor  Is 
critical  on  prompt  neutons  alone.  The  simulated  transients  Imposed  on 
the  model  were  limited  to  ensure  that  excessive  amounts  of  reactivity 
were  not  Introduced. 

B.  REACTIVITY  FEEDBACK  MECHANISMS 

In  the  zero  power  reactor  point  reactor,  the  power  level  Is  assumed 
to  be  so  low  that  It  does  not  affect  the  reactivity,  thus  there  are  no 
feedback  effects.  However,  for  a  reactor  operating  at  a  useful  power 
level,  feedback  effects  do  exist  and  the  reactivity  becomes  an  Implicit 
function  of  the  reactor  power  level  (or  neutron  flux).  This  dependence 
arises  since  reactivity  depends  on  macroscopic  cross  sections  which  In¬ 
volve  the  atomic  number  densities  of  material  In  the  core: 

2  -  No¬ 
where  £  *  Macroscopic  cross  section  (cm-1) 

N  ■  Atomic  number  density  (atoms/cm5  ) 

<T  *  Microscopic  cross  section  (cm  *  ) 

The  atomic  number  density  can  depend  upon  the  reactor  power  level 
since  the  concentrations  of  certain  nuclei  are  constantly  changing  due 
to  neutron  Interactions.  Material  densities  also  depend  upon  tempera¬ 
ture  which  Is  a  function  of  reactor  power  level  and  hence  the  fluxJ 

^uderstadt,  J.  J.  and  Hamilton,  L.  J.,  Nuclear  Reactor  Analysis, 

1st  ed.,  Wiley,  1976. 


The  two  feedback  mechanisms  considered  In  the  model  were  Xenon-let  and 
moderator  temperature. 

Reactivity  can  also  be  changed  directly  by  an  external  source  such 
as  control  rods  containing  a  neutron  absorbing  material.  Thus  the  over¬ 
all  reactivity  can  be  written  as: 

*  tyr  (12) 

where  A«  Fx  P 
&Tm  Fj  P 

■  Overall  core  reactivity 
•  Externally  added  reactivity 
4>*  »  Xenon-135  feedback  reactivity 
&T  ■  Moderator  temperature  feedback  reactivity 
Fx  -  Xeonon  transfer  function 
Fy  ■  Moderator  temperature  transfer  function 

Equation  (12)  with  the  previously  developed  zero  power  reactor 
transfer  function,  given  by  equation  (11),  are  incorporated  In  the 

0 

block  diagram  below: 


Figure  3. 


Block  Olagram  of  the  Reactor  Transfer  Function 
with  Xenon-135  and  Moderator  Temperature 
Feedback  Loops,  and  External  Reactivity 
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The  next  step  In  the  model  development  was  to  derive  the  Xenon-135 


and  moderator  temperature  transfer  functions. 

1 .  Xenon-135  Transfer  Function 

Xenon-135  Is  the  most  significant  fission  product  poldon  because 
of  Its  enormous  thermal  neutron  absorption  cross  section  and  relatively 
large  fission  yield.  This  Isotope  Is  not  only  produced  directly  from 
fission  but  also  from  the  decay  of  other  fission  products  as  shown  In 
Figure  4. 


Figure  4.  Xe-135  Decay  Scheme 


Since  the  /3“decay  of  Iodine-135  (135I)  and  Xenon-135  (135Xe),  with 

the  largest  half-lives,  are  the  controlling  steps  In  this  decay  scheme. 

It  was  simplified  by  making  the  following  assumptions: 

1 35 

1)  All  I  Is  produced  directly  from  fission  (the  production  of 

Antimony-135  (133Sb)  and  subsequent  decay  to  ^35I  is  considered 
instantaneous). 

2)  The  short  lived  metastable  135mXe  Is  Ignored. 

1 35 

3)  The  removal  of  I  by  neutron  absorption  Is  negligible  for 

the  neuton  flux  levels  used  In  the  model  (10^  neutrons/cm2sec). 


4)  All  the  135 1  decays  to  1j3Xe. 


With  these  assumptions  the  effective  decay  scheme  is  shown  In  Figure 


Figure  5.  Simplified  Xe-135  Decay  Scheme 

135 

Using  this  decay  scheme,  the  resulting  rate  equations  for  I  and 
Xe  are: 

tfxSfdCO  -Axl<t>  (13) 

Yx^eJ(t)  +XiIte)-Axx<t)-<Jax^Ct)X(t)  (14) 

1 35  3 

where  I  *  I  number  density  (atoms/cm  ) 

X  *  135Xe  number  density  (atoms/cm3) 

Zf  *  Macroscopic  fission  cross  section  (cm-1) 

2 

0  *  Integrated  one  group  neutron  flux  (neutrons/cm  sec) 
Ifj  *  Effective  ^3®I  fission  yield 

■  Effective  ^Xe  fission  yield 

Aj  ■  decay  constant  (sec  ) 

1 35  -1 

Ax  *  Xe  decay  constant  (sec  ) 

■  ^3^Xe  microscopic  thermal  neutron  absorption  cross 
a  section  (cm2) 


Again,  using  the  first  order  perturbation  technique,  let 
I(t)  -  io  +  <ft(t) 

X(t)  *  XQ  +  <A(t)  (15) 

0(t)  -  0Q  +  oty(t) 

where  the  zero  subscript  denotes  an  initial  steady  state  value  and 
the  delta  prefix  denotes  a  small  perturbation  about  this  value. 

Upon  substituting  equations  (15)  into  equations  (13)  and  (14),  and 

taking  the  Laplace  transform,  the  relationship  between  the  perturbation 
1 35 

in  Xe  and  0  is  derived  (See  Appendix  A): 


_  afxSf-csXXols+XtQfrEf  X.) 

(o5*d«  +  *K*>Ax>st  *r(<Ja*  to*  A*) 


=  6jcts) 


(16) 
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The  change  in  reactivity  caused  by  a  small  perturbation  in  Xe 
concentration  is  also  derived  in  Appendix  A: 


where 


_  1 

"  Xo  +  *Lu 
<5* 


1 35 

X  *  Equilibrium  Xe  number  density  before  the 
perturbation  (atoms/cm3) 

235  3 

U  3  Uranium-235  (  U)  number  density  (atoms/cm  ) 

3  Microscopic  thermal  neutron  absorption  cross 
section  of  2^ U  (cm2) 


(17) 


A  perturbation  in  neuton  flux  (  c^0)  is  directly  proportional  to  a 
perturbation  in  reactor  power  (  ^*P)  as  shown  by: 
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where  *  Energy  released  per  fission 

*  Macroscopic  fission  cross  section  of  U  (cm  ) 


or 

</^(s)  _  1 

cZ'P(s)  Ef  Sf  =  K 

The  product  of  equations  (17),  (16),  and  (18)  yields: 


(18) 


ofocfrl  e^XQl  </)0(s)  _  <7^<(s) 
</*X(s)  ^(s)  (/‘PCS)  c/’P(s) 


«<XGX<S)K  =  Fx(s) 


(19) 


Equation  (19)  is  the  transfer  function  relating  reactor  power  and 
135 

Xe  feedback  reactivity.  This  equation  is  shown  in  block  diagram 
form  in  Figure  6  where  it  has  been  incorporated  with  the  previously  de¬ 
rived  zero  power  transfer  function  given  by  equation  (11). 


Figure  6.  Block  Diagram  of  the  Reactor  Transfer 
Function  with  Xenon-135  Feedback  and 
External  Reactivity 


2.  Moderator  Temperature  Feedback 

Reactivity  feedback  from  changes  In  moderator  temperature  occurs 
as  a  result  of  changes  with  moderator  density.  The  density  is  also  a 


24 


function  of  pressure,  however,  the  pressure  coefficient  of  reactivity  is 
typically  two  orders  of  magnitude  smaller  than  the  temperature  coeffi¬ 
cient,  and  was  therefore  not  considered  in  the  model  development.  The 
moderator  density  affects  the  moderator  number  density  and  hence  the 
macroscopic  scattering  cross  section  of  the  moderator. 

The  primary  mechanism  for  the  thermal ization  of  the  prompt  and  de¬ 
layed  neutrons  is  by  elastic  scattering  interactions  with  the  moderator 
nuclei.  Hence,  variations  in  the  macroscopic  scattering  cross  section 
will  affect  the  rate  at  which  neutrons  become  thermalized.  Changes  in 
the  thermal ization  rate  affect  the  fission  rate,  or  reactor  power  level. 
Power  level  changes  affect  the  moderator  temperature,  thus  a  feedback 
loop  is  created. 

Because  the  moderator  is  also  the  coolant  in  the  model,  its  temper¬ 
ature  is  not  only  a  function  of  reactor  power  but  also  a  function  of  the 
heat  transfer  process  occurring  in  the  heat  exchanger,  thus  complicating 
the  feedback  loop. 

Therefore,  in  order  to  develop  this  feedback  mechanism  analytically, 
it  is  first  necessary  to  model  the  plant's  heat  transfer  processes. 

This  thermal  analysis  is  done  in  the  following  section. 

Another  temperature  feedback  mechanism  is  the  broadening  of  the 
Uranium-238  resonance  absorption  cross  section  for  neutrons  with  increas¬ 
ing  temperature.  Because  a  highly  enriched  core  was  assumed  in  the 
model,  with  little  Uranium-238,  this  effect  was  not  considered. 
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C.  THERMAL  ANALYSIS 


1 .  Reactor  Heat  Transfer  Function 

A  lumped  parameter  model  was  assumed.  This  simplification  pro¬ 
vided  a  set  of  ordinary  differential  equation  which  were  sufficiently 
accurate  for  the  simulated  normal  operating  transients.  In  the  lumped 
parameter  model,  heat  transfer  in  the  reactor  was  assumed  to  occur  at  a 
single  point.  Thus,  spatial  variations  were  neglected.  The  core  was 
considered  to  be  a  homogenized  mixture  of  the  uranium  alloy  fuel,  the 
fuel  cladding,  and  other  structural  materials,  with  a  constant  thermal 
capacity. 

The  equation  for  the  heat  flow  from  the  core  was  obtained  from  a 
basic  heat  balance.  The  heat  generated  from  fission  equals  the  heat 
required  to  change  the  temperature  of  the  core  materials  plus  the  heat 
transferred  to  the  coolant.  On  a  per  unit  time  basis,  the  heat  balance 
is: 

dT-(t)  r  “I 

p(t)  >  CF  [  TpCt)-  TAV(t)  J  (20) 

where  P(t)  ■  Total  Power  generated  in  the  core  (Btu/sec) 

Tp(t)=  Average  temperature  of  the  core  materials  (°F) 

Tfty(t)  *  Average  reactor  coolant  temperature  (°F) 

Cp  »  Total  thermal  capacity  of  the  core  materials  (Btu/°F) 
hpm  *  Total  heat  transfer  coefficient  (Btu/°F  sec) 

A  similar  heat  transfer  balance  must  also  hold  for  the  heat  being 
transferred  to  the  coolant  and  transported  out  of  the  core.  This  heat 
is  the  last  term  of  equation  (20)  and  is  transferred  to  the  coolant  as: 
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a  Cm  si3^(fe>  +  wmC  [  TH#Co (21) 

where  CM  *  Total  thermal  capacity  of  coolant  in  the  core 

M  (Btu/°F) 

m  =  Coolant  mass  flow  rate  (Ibm/sec) 
m 

Tu  *  Average  reactor  coolant  outlet  temperature  (°F) 

Ho 

T_  *  Average  reactor  coolant  inlet  temperature  (°F) 
ci 

C  =  Specific  heat  of  reactor  coolant  (Btu/lbm°F) 

For  simplification,  the  average  reactor  coolant  temperature  is 
assumed  to  be  given  by: 


Tav  (t)  =  £*!«*(*)  +  Th0  (*>]  /  2 
Rearranging  equation  (20)  yields: 


where 


Tp(t)  i-  rt 


dT|»ttl  _ 
d't 


TavU)  + 


J1 

CF 


PCt) 


(22) 


(23) 


Substituting  equation  (22)  into  equation  (21)  to  eliminate  T^yields: 


TawM  ♦  4  -  tku  +  *J|tsw 


where 


To  *  StL.  t  and 

m„C  * 


(24) 
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P(t)  3  PQ  +  <^P(t) 

Tp(t  )=  Tp  +  ^Tp(t) 

0 

TAV(t)*  TAV  +  (25) 

o 

V’-V  Ate'1' 

0  0 

Tc<*)*Tci*  A<*> 

1  0  0 

where  the  zero  subscript  denotes  a  steady  state  value  and  the 
prefix  a  small  perturbation  about  this  value. 

Substituting  equations  (25)  Into  equations  (22),  (23),  and  (24), 
eliminating  the  average  temperatures  Tp(t)  and  T.y(t) ,  and  taking  Laplace 
transforms,  the  following  relationship  between  reactor  coolant  inlet  and 
outlet  temperature  perturbations  is  derived  (see  Appendix  B): 


”  {  ^  S*  *  [  If  ]  S  "  1  }  t  *0^(5) 

^Th<i(sJ  =  - - - 

where 

*  iik  =  ,x  - 

Cp7* 


(26) 


This  equation  is  the  reactor  heat  transfer  function  for  the  reactor 
coolant  outlet  temperature  as  a  function  of  the  reactor  coolant  inlet 
temperature  and  reactor  power. 

A  block  diagram  representation  of  this  transfer  function  is  shown 


in  Figure  7. 


y>W 


Figure  7.  Block  Diagram  of  the  Reactor  Heat 
Transfer  Function 


In  Figure  7 

GC(S)  - 

G„(S)  ■ 


2.  Heat  Exchanger  Transfer  Function 

As  in  the  derivation  of  the  reactor  heat  transfer  function,  a 
lumped  parameter  model  was  assumed  for  the  heat  exchanger.  Two  points 
of  energy  storage  were  assumed,  the  primary  coolant  water  and  the  water 
on  the  secondary  side  of  the  heat  exchanger.  The  thermal  capacity  of 
the  heat  exchanger  metal  was  included  with  that  for  the  secondary  water, 
since  the  thermal  resistance  on  the  primary  side  Is  predominant.  A 
further  assumption  was  made  that  the  time  spent  by  the  primary  coolant 
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while  It  passed  through  the  heat  exchanger  was  negligible  In  comparison 
with  the  time  spent  In  the  primary  piping. 

The  state  of  the  steam  produced  on  the  secondary  side  of  the  heat 
exchanger  was  assumed  to  always  be  dry  and  saturated  and  that  the 
secondary  water  Is  always  at  the  saturation  temperature  for  the  existing 
pressure.  These  assumptions  were  justified  because  moisture  separators 
and  recirculation  can  provide  high  quality  steam  and  preheating  of  the 
feedwater. 

With  these  assumptions,  the  following  equations  were  obtained  from 
a  heat  balance  per  unit  time: 

-Tc*  (fc)  +  H-n*  [l*ve(0  -  TsU:)^  (27) 


£lAV8(tl  —  TjW  ^  —  C§  +  ft  Ot) 


(28) 


where 


THi(t) 

TC0^> 

TAVB(t) 

Ts(t) 


’TM 

Si 

Ce 


PL(t) 


m_ 


Heat  exchanger  coolant  inlet  temperature  (°F) 
Heat  exchanger  coolant  outlet  temperature  (°F) 
Heat  exchanger  average  coolant  temperature  (°F) 

Saturated  steam  temperature  (°F) 

Total  heat  transfer  coefficient  (BTU/°F.sec) 

Total  thermal  capacity  of  coolant  In  the  heat 
exchanger  (BTU/°F) 

Total  thermal  capacity  of  heat  exchanger  metal 
and  secondary  water  and  steam  (Btu/°F) 

Power  delivered  by  the  heat  exchanger  (BTU/sec) 
Coolant  mass  flow  rate  (Ibm/sec) 
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The  total  heat  transfer  coefficient,  a  function  of  the  primary  coolant 
flow  rate  (a  constant  In  the  model)  and  the  heat  transfer  characteristics 
of  the  heat  exchanger,  was  assumed  to  be  constant.  Also,  the  time  delay 
In  transferring  heat  across  the  heat  exchanger  tubes  was  neglected. 

This  changed  the  shape  of  the  Initial  thermal  transient  but  had  little 
effect  on  the  basic  dynamics  of  the  secondary  loop. 

Again,  as  a  simplification,  the  average  temperature  of  the  coolant 
In  the  heat  exchanger  was  assumed  to  be  given  by 

S  ■  1  <th,  *  v  (29) 

The  power  delivered  by  the  steam  generator  is  proportional  to  the 
product  of  the  steam  flow  rate  and  the  difference  in  enthalpy  between 
the  steam  and  feedwater 

where  m$  *  Steam  flow  rate  (Ibm/sec) 

Hs  *  Saturated  steam  enthalpy  (Btu/lbm) 

Hp  «  Feedwater  enthalpy  (Btu/lbm) 

This  equation  assumes  the  steam  and  feedwater  flow  rates  are  always 
equal  and  thus  neglects  any  Instabilities  In  the  secondary  steam. 

As  discussed  In  Reference  4,  the  enthalpy  of  saturated  steam  Is 
nearly  constant  over  a  wide  range  of  pressure,  varying  from  1198  Btu/lbm 
to  1204  Btu/lbm  over  a  pressure  range  from  200  to  800  pslg.  The  feed- 
water  enthalpy .which  depends  on  condenser  pressure.  Is  usually  between 
50  and  100  Btu/lbm.  Thus  the  enthalpy  difference  may  be  regarded  as  a 


constant,  and  the  power  delivered  by  the  heat  exchanger  Is  directly 
proportional  to  the  steam  flow  rate. 

The  impedance  to  steam  flow  caused  by  the  turbine  Is  nearly  Indepen 
dent  of  turbine  speed.  If  constant  backpressure  Is  assumed,  the  steam 
flow  rate  Is  directly  proportional  to  the  throttle  opening  at  a  given 
pressure.  Thus, 

«s  *  y*sA 

2 

where  «  Saturated  steam  pressure  (Ibf/in  ) 

A  *  Proportionality  factor  which  is  a  function 
of  the  throttle  setting  TJm1n£ 

1 57  sec 


The  numerous  assumptions  made  in  the  heat  exchanger  model  development 
limit  the  accuracy  of  the  resulting  equations.  However,  the  errors  in¬ 
volved  in  these  assumptions  are  usually  less  than  the  amount  of  uncertainty 
in  the  engineering  valves  of  the  coefficients  used  in  the  equations. 

Recalling  the  previous  assumption  that 


dTAVB(0 

dt 


=  0  ,  and 


Substituting  equation  (29)  into  equation  (27),  yields: 


'"mC  *  ***"{? [ Thi <*>■•■  Tco(t)]  -  Ts<t)| 

Solving  equation  (30)  for  Tr  ,  gives: 

S 

•T  ,  2Ts(b)  -  Th*  (1-Ki) 

lcev«l  =  - - - 

ri*Ki) 


[30) 


(31) 
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where 


z 


K 


1 


Ktm 


Substituting  equation  (29)  into  equation  (28)  gives: 

Kr*  (t  [t*HW  +  TSiW  ]  -  Ts<0 1  *  Cs  +  PUt) 

Rearranging  this  expression. 


Cs  cjr8 
Ktm  dt 


+  Ts(e)  —  +  Tc«<o]  +  R-CfcJ 


(32) 


Let  Thi  (t )  =  Th  +  (/^(t) 

Tc  (t)  =  Tc  +  ^Ts(t)  (33) 

<3  °o 

Ts(t)  *  Tx0  +  ^z{t) 

pL(t)  -  PL  +  </*pL(t) 

where  the  lowest  zero  subscript  denotes  a  steady  state  value  and 
the  delta  prefix  a  small  perturbation  about  this  value. 

By  substituting  equations  (33)  into  equations  (31)  and  (32),  and 
taking  the  Laplace  transform,  the  following  set  of  equations  is  derived 
(see  Appendix  C): 

t  -  fa  A(s) 

=  - -i - * -  (34) 

Ts  +  1 


y’TcoCS)  « 


g^*ls(5)  -  [l-  *il 

1+Ki 


where 


(35) 


i 
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These  equations  are  represented  In  a  block  diagram  In  Figure  1. 


Figure  8.  Block  Diagram  of  the  Heat  Exchanger 
Transfer  Function 


With  ^Tu  and 
H1 

algebraic  loop  In 


^PL  as  Inputs,  equations  (34)  and  (35)  form  an 
^Tj.  and  These  two  equations  were  solved 


with  the  use  of  an  Inherent  functional  routine  available  In  the  CSMP-III 


language  which  was  used  to  formulate  the  model. 
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3.  Primary  Piping  Transfer  Functions 

While  circulating  through  the  primary  loop,  the  coolant  under¬ 
goes  mixing  and  transport  delay  effects.  For  example,  a  transient  in 
the  coolant  temperature  at  the  heat  exchanger  outlet  does  not  appear  at 
the  reactor  Inlet  until  sometime  later.  Where  there  are  volume  or  flow 
direction  changes,  as  In  the  reactor  coolant  Inlet  plenum,  mixing  occurs 
causing  a  smearing  of  a  temperature  transient.  Both  of  these  effects 
were  approximated  in  the  model  by  combinations  of  two  types  of  time 
delays:  a  pure  transport  delay  and  a  simple  mixing  delay. 

Assuming  no  mixing  in  the  primary  piping  and  no  heat  loss  with 
perfectly  insulated  pipes,  pure  transport  delays  are  encountered  In  the 
piping  runs  between  the  reactor  and  heat  exchanger.  With  these  assump¬ 
tions,  the  temperatures  Involved  in  the  transfer  of  heat  between  the  re¬ 
actor  and  heat  exchanger  can  be  given  as 

y*>  *th0  "-v 

where  the  inlet  coolant  temperature  to  the  heat  exchanger  has  the 

same  form  as  the  outlet  temperature  of  the  reactor  TH  after  a  fixed 

o 

transport  delay  Tj.  This  transport  delay  can  be  approximated  by  the 

following  differential  equation  which  Is  derived  in  Appendix  D.  i 

T„  *  I',  ^1  ’  T  (36) 

i  at  o 

Similarly,  the  transport  delay  from  the  heat  exchanger  outlet  to  j 

the  reactor  Inlet  plenum  Is:  j 

Te1 +  ^  dTcip(t)  »  T-  (t)  (37) 

1p  St —  c° 
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where 


Tr  (t)  *  Reactor  Inlet  plenum  coolant  temperature  (°F) 
Ip 

Tr  (t)  *  Heat  exchanger  outlet  coolant  temperature  (°F) 

Lo 

This  can  be  approximated  by 

Tc^(fc)  +  h  ^  *  rc°(t)  {3?) 


As  in  the  previous  development,  let 


T 

T 

T 

T 


(t) 


ip 


(t) 


(t) 


'ip 


(t) 


<Ah  (t) 

<fTu  (t) 

Hip 

(O 

Co 

<AC  (t) 

Cip 


(38) 


where  again  the  lowest  zero  subscript  denotes  a  steady  state  value 
and  the  delta  prefix  a  small  perturbation  about  this  value. 

After  substituting  equations  (38)  into  equations  (36)  and  (37), 
taking  the  Laplace  transform,  the  following  transfer  functions  are 
derived  (see  in  Appendix  0). 


</*Th^(5)  _  1 

<^T»loCS>  1+^3 s 


(39) 


/Tc^Cs)  _  1 

<y*7c0CS)  1  NS 
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Combined  mixing  and  transport  effects  are  encountered  at  the  reactor 
and  heat  exchanger  inlet  and  outlet  coolant  plenums.  For  simplification, 
only  the  mixing  effects  at  the  inlet  plenums  were  considered  in  the 
model.  Assuming  perfect  mixing  and  after  performing  a  heat  balance  on 
the  plenum  concerned,  the  combined  mixing  and  transport  effects  are  ex¬ 
pressed  by  differential  equations  developed  in  Appendix  D  of  the  form 

♦  To(t>  -  T1(t) 

where 

T\(t)  *  Plenum  inlet  temperature  (°F) 

T0(t)  *  Plenum  outlet  temperature  (°F) 

*f  *  Mixing  delay 

As  before,  by  using  a  small  perturbation  technique,  and  taking  the 
Laplace  transform,  the  following  transfer  function  is  derived  (see 
Appendix  D): 


=  — V-  (41) 

</*7i^(s)  1  ns 

Block  diagram  representations  of  the  coolant  piping  transport  and 
mixing  transfer  functions  is  shown  in  Figure  9. 

D.  OVERALL  PLANT  BLOCK  DIAGRAM 

The  transfer  functions  developed  previously  for  the  zero  power  point 
1 35 

reactor,  Xe  feedback,  reactor  and  heat  exchanger  heat  transfer,  and 
the  primary  piping  were  Interconnected  resulting  in  the  overall  model 
block  diagram  shown  in  Figure  10.  The  moderator  temperature  feedback 
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</>(*) 


Figure  9.  Block  Diagram  of  the  Transport 
and  Mixing  Delay  Functions 


loop  is  seen  to  consist  of  the  heat  transfer  process  and  primary  piping 

transfer  functions  which  yield  o^T-  and  .  These  temperatures 

Ci  "o 

are  summed,  then  halved,  yielding  T which  is  then  multiplied  by 
the  negative  temperature  coefficient  °*T,  generating  the  moderator 
temperature  feedback  reactivity  fi  j. 

E.  MODEL  CONSTANTS 

The  following  plant  parameters  from  the  Shippingport  Atomic  Power 
Station  were  obtained  from  Reference  6  and  used  in  the  model  development. 
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1 .  General  Parameters 

a.  Reactor  thermal  power 

b.  Reactor  coolant  system  pressure 

c.  Reactor  coolant  average  temperature 

d.  Steam  pressure  at  full  load 


231  MW 
2000  psla 
523°  F 
600  psla 


2.  Reactor  Coolant  System 

a.  Reactor  coolant  flow  rate 

b.  Reactor  coolant  outlet  temperature 

c.  Reactor  coolant  inlet  temperature 

d.  Coolant  volume  in  core 


6280  Ibm/sec 
538°  F 
508°  F 
103  ft3 


3.  Reactor  Core 


a.  Configuration 

Right  cylinder 

b.  Size 

6.8  ft.  dia  X  6  ft.  high 

c.  Fuel  load  U235  (seed) 

75  kg 

d.  Composition  (seed) 

1 )  Water 

43.5  v/o 

2)  Fuel  alloy 

30.0  v/o 

3)  Zircalloy 

34.1  v/o 

e.  Control  Rods 

1 )  Total  rod  worth 

0.256  k 

2)  Scram  time 

0.35  sec  time  delay 

1 .0  sec  rod  drop 

Nuclear  Data 

2  X  1014  n/cm2  sec 
5.6  X  10-5  sec 


a.  Thermal  neutron  flux 

b.  Prompt  neutron  lifetime 
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c.  Effective  delayed  neutron  fraction  0.0077 

d.  Temperature  coefficient  of  reactivity  -3.1  X  10”4  ^*k/ °F 

5.  Reactor  Protection  Setpoints 

a .  Scram 

1)  High  reactor  power  138% 

2)  High  reactor  coolant  outlet  temperature  550°F 

b.  Cutback 

1)  High  reactor  power  114% 

2)  High  Startup  rate  1.74  Decades  per  min. 
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IV.  RESULTS 


A.  INTERACTIVE  PROGRAM 

Figures  10-16  are  examples  of  the  Interactive  program  output.  As 
shown  in  Figure  10,  the  user  is  first  given  a  description  of  the  purpose 
of  the  overall  program  and  general  instructions  for  input  entry.  As 
shown  in  Figure  11,  the  user  is  then  prompted  to  enter  an  initial  steady 
state  power  level,  after  which  the  corresponding  major  plant  parameters 
are  displayed.  The  user  is  then  prompted  to  choose  the  type  of  simula¬ 
tion,  either  plant  transient  at  power,  or  post-shutdown  Xenon-135 
behavior.  The  user  must  be  at  the  Tektronix  4012  console  in  the  computer 
center  if  the  latter  is  chosen. 

If  the  post-shutdown  Xenon-135  behavior  is  chosen  to  be  examined, 
graphs  such  as  those  shown  in  Figure  13  are  generated  and  displayed  on 
the  Tektronix1 s  screen.  In  addition  to  showing  the  time  response  of  the 
Xenon,  the  time  of  its  peak  and  associated  maximum  reactivity  are  dis¬ 
played. 

If  the  plant  transient  simulation  is  chosen,  the  user  is  then  prompted 
to  choose  the  mechanism  for  initiating  the  transient,  either  control  rod 
movement  or  turbine  load  change.  In  Figure  14,  the  user  has  chosen  the 
former.  The  program  then  prompts  the  user  to  enter  the  direction  and 
time  length  of  the  control  rod  movement.  In  Figure  15,  the  turbine  load 
change  was  chosen  and  the  program  request  the  final  turbine  load.  In 
either  case,  a  summary  of  the  transient  Inputs  is  displayed  with  a 
notice  that  the  Interactive  portion  is  complete. 
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Figure  11.  Interactive  Program,  Introduction  and  Instructions 


ENTER  AN  INITIAL  STEADY  STATE  POWER  LEVEL  BETWEEN  10.  AND  100.  PERCENT 
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ENTER  THE  INITIAL  REACTIVITY  CHANGE  MECHANISM 

1.  BANK  CONTROL  ROD  CHANGE 

2.  TURBINE  LOAD  CHANGE 
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Figure  14.  Interactive  Program,  Control  Rod  Movement  Inputs 
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Figure  15.  Interactive  Program,  Turbine  Load  Change  Inputs 
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Figure  16  Is  an  example  of  the  Interactive  program's  logic  routines 
which  will  only  permit  inputs  within  the  requested  ranges  to  be  accepted 
for  the  simulation. 


8.  SIMULATION  PROGRAM 


The  program  generates  Versatec  plots  showing  the  plant's  power,  re¬ 
activity,  and  temperature  transient  responses  for  the  user's  Inputs. 
Figures  17-20  are  composites  of  these  plots. 

In  Figure  17,  the  transient  was  Initiated  by  a  simulated  20  to  40 
percent  ramp  change  in  turbine  load  at  1/2  percent  per  second.  In 
Figure  18,  the  transient  was  initiated  by  a  simulated  60  to  40  percent 
change  in  turbine  load  in  the  same  manner.  Both  of  these  figures  show 
the  characteristic  power  demand  following  and  inherent  stability  response 
typical  of  a  PWR  plant  with  a  constant  average  temperature  program  and 
negative  temperature  coefficient. 

The  inherent  stability  feature  is  further  displayed  in  Figures  19 
and  20.  In  Figure  19,  the  transient  is  initiated  by  a  simulated  10  second 
Inward  movement  of  the  control  rods  at  a  reactivity  insertion  rate  of 
1.25  x  10  d k  per  second  from  an  initial  50  percent  power  level.  In 
Figure  20,  the  transient  response  to  an  outward  control  rod  movement 
with  the  same  parameters  is  shown.  As  seen  in  both  cases,  the  reactor 
power  stabilizes  about  the  initial  turbine  load  and  the  average  reactor 
coolant  temperature  stabilizes  at  a  new  level  to  compensate,  via  the 
negative  temperature  coefficient,  for  the  reactivity  inserted  by  the 
control  rods. 
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Figure  16.  Interactive  Program,  Logic  Routine  Example 
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Figure  18.  Power,  Reactivity,  and  Temperature  Response 
to  a  60  to  40%  Turbine  Load  Change 


Figure  19.  Power,  Reactivity,  and  Temperature  Response 
to  a  10  second  Inward  Control  Rod  Movement 
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SIMULATION  OF  A  PRESSURIZED  WATER  REACTOR 
POWER  AND  LOAD  RESPONSE  legend 
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Figure  20.  Power,  Reactivity,  and  Temperature  Response 
to  a  10  second  Outward  Control  Rod  Movement 


V.  CONCLUSIONS  AND  RECOMMENDATIONS 


The  model  responds  to  normal  operating  transients  with  the  charac¬ 
teristics  of  a  PWR  plant.  The  interactive  program  provides  the  user  with 
the  facility  to  initiate  the  simulation  or  examine  a  real  time  display 
of  post-shutdown  Xenon-135  buildup  and  decay. 

As  previously  discussed,  the  purpose  of  this  work  was  to  develop  a 
learning  device  for  nuclear  engineering  students  at  the  Naval  Postgrad¬ 
uate  School.  The  result  was  a  relatively  simple  model  of  a  PWR  power 
plant  and  a  simulation  program  that  suffers  from  a  long  execution  time. 

The  model's  simplicity  resulted  mainly  from  the  consideration  of 
an  effective  single  group  of  delayed  neutron  precursors  in  the  reactor 
kinetic  equations  and  the  use  of  lumped  parameter  analysis  in  the  plant's 
heat  transfer  processes.  The  model's  sophistication  could  be  increased 
by: 

1.  Expanding  the  reactor  kinetic  equations  to  consider  six 
delayed  neutron  precursor  groups. 

2.  Developing  a  multi-section  heat  transfer  model  for  the 
reactor  and  the  heat  exchanger. 

3.  Incorporating  variable  reactor  coolant  mass  flow  rate  and 
heat  transfer  coefficients. 

These  additional  features  will  also  increase  the  complexity  of  the 
simulation  program,  thus  aggravating  the  long  run  time  problem.  This 
concern  might  be  eased  by  taking  a  different  approach  to  the  model's 
formulation  than  the  transfer  function  method.  While  this  method  was 
readily  coded  using  the  CSMP-III  language,  physical  insight  to  the  plant 
dynamic  processes  was  lost  when  the  Laplace  transform  and  subsequent 
grouping  of  constant  terms  was  performed.  The  use  of  state  variable 
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theory  might  not  only  allow  this  physical  appreciation  to  be  retained 

but  also  result  in  shorter  program  run  times. 

Regardless  of  the  formulation  method,  a  real  time  response  to  the 

simulation  can  probably  only  be  achieved  by  using  an  analog  computer. 

The  wide  range  of  time  constants  associated  with  the  equations  (prompt 

-5 

neutron  life-times  on  the  order  of  10  seconds  to  Xenon  decay  half 
lives  on  the  order  of  hours)  require  a  small  numerical  integration  time 
interval  over  a  large  time  period.  The  result  is  a  prohibitively  long 
run  time  on  a  digital  computer  for  an  interactive  program. 

It  is  hoped  that  refinement  of  the  model  will  be  continued.  While 
the  existing  model  does  reflect  the  general  characteristics  of  a  PWR 
plant  to  normal  operating  transients,  it  does  exhibit  relatively  large 
power  over/under  shoots  in  response  to  turbine  load  changes  as  seen  in 
Figures  17  and  18. 
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APPENDIX  A 


DEVELOPMENT  OF  THE  XENON  FEEDBACK  TRANSFER  FUNCTION  COMPONENTS 


1 .  Derivation  of  Gx (S ) 

The  previously  developed  rate  equations  for  135I  and  135Xe  from  the 
1 35 

effective  Xe  decay  scheme  shown  in  Figure  5  are: 


dl 


-  &  Ilf  6  ct)  -  Axr ct) 


dx 


=  *x  2f  6(t)  +  Ail  (t)  -  AxX(^l  —  &(t)  X(t) 


135  3 

where  I(t)  =  I  number  density  (atoms/cm  ) 

X(t)  =  133Xe  number  density  (atoms/cm3) 


(A  1) 
(A  2) 


Sf  *  Macroscopic  fission  cross  section  of  (cm-1) 

0  a  Average  integrated  one  group  flux  (neutrons/cm  sec) 
^1  =  Effective  135I  fission  yield 

^X  *  Effective1 35Xe  fission  yield 
^1  *  135I  decay  constant  (sec-1) 

^X  =  135Xe  decay  constant  (sec-1) 


-1 


1 35  2 

X  *  Xe  thermal  neutron  absorption  cross  section  (cm) 


<S 


Let  I(t)  =  IQ  +  </l(t) 
X(t)  =  XQ  +  d\{ t) 
0(t)  =  0Q  +  </0{t) 


(A  3) 


where  the  zero  subscript  denotes  an  initial  steady  state  value  and 
the  delta  prefix  denotes  a  small  perturbation  about  this  value. 
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With  the  reactor  at  an  initial  steady  state  with  a  flux  level  0  ,  the 
1 35  135 

equilibrium  values  of  I  and  Xe  are  found  by  setting  their  time 

dependence  in  equations  (A  1)  and  (A  2)  equal  to  zero.  Thus, 


,  .  VxZfzSo 
A, 

^  Cta  +  Ax 

*X  +  <3"d*  do 

Substituting  equation  (A  4)  into  equation  (A  5) 


(A  4) 

(A  5) 


u»+yx)  Zt<£ o 

Ax  +  O5x^o 


(A  6) 


Substituting  equations  (A3),  (A  4),  and  (A  6)  into  equation  (A  1) 
and  (A  2),  neglecting  the  <^0  cA(  term  (first  order  approximation) 


=  J/rSfc/fc-Arcrtr  (A  7) 

jL  »Ax  =  AzSl  t  ( *x  Zp  -Oj'x.)  c/tf-CA*  +03*4.) cfx  (A  8) 

Taking  the  Laplace  transform  of  equations  (A  7)  and  (A  8)  and  solving 
for  c/l(S)  and  cfx(S) 


cPl(s)  - 


Vr  Spc/^Cs) 

sTai 


c/XCS)  = 


Ar<^I  (s)  <■  C  2P  -  Qa*x<>) 
(st  AX  +aax&,) 


(A  9) 

(A10) 


J 
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Substituting  equation  (A  9)  into  equation  (A10) 


,  .  -csaxX«)o^CsO 

oXOi  *  " 

(SiKx)  (  St\nt  Oa^d) 

<^x(sl 

Expanding,  collecting  terms,  and  solving  for 


JXCS)  ^  (J±Zf  -<3a*Xe)s+  AxC^rSp  -Q5X^»)  __  G^(s) 

^ ®  s*  +  (0£X4o  +  Ax  ♦■Ax)  s  +  (<£*<*«  +  Ax) 

This  is  equation  (16)  on  page  23. 


2.  Derivation  of  °*x 

The  effective  core  neutron  multiplication  factor  k  is  given  by 
the  familiar  "six  factor  formula"  found  in  the  literature 


k  -  pna  "tnl  '*>’> 

where  ^  *  Thermal  fission  factor  of  the  fule 
6  *  Fast  fission  factor  of  the  fuel 

p  =  Resonance  escape  probability 
f  =  Thermal  utilization  factor  in  the  core 
PFNL  *  Fast  neutron  non-leakage  probability 
PTNL  *  Thermal  neutron  non-leakage  probability 


The  core  reactivity  is  defined  as 


(A  12) 


Consider  a  reactor  at  an  initial  critical  steady  state  condition 
1 35 

with  an  equilibrium  Xe  concentration  Xo.  By  the  definition  of 
criticality 
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In  this  condition 


ko  3  1  {A  13) 

/Oq  *  o 


(A  14) 


where 


Fuel  macroscopic  thermal  neutron  absorption  cross 
section 


rcore 


*  Initial  core  macroscopic  thermal  neutron  absorption 
cross  section 


2 1°  - 


135 

Initial  Xe  macroscopic  thermal  neutron  absorption 
cross  section 


and  the  macroscopic  thermal  neutron  absorption  cross  sections  of  the 
moderator  and  other  core  materials  has  been  neglected,  a  reasonable  ap¬ 
proximation  in  a  highly  enriched  core. 

1 35 

Now  impose  a  small  perturbation  in  the  Xe  concentration  on  the 

1 35 

core.  Neglecting  the  effect  of  the  Xe  perturbation  on  neutron  leakage, 
a  reasonable  assumption  for  a  large  reactor, 

X  *  Xo  +  </*X 
k  ■  ko  +  cfk 
f  a  fO  +  cTf 

/>  = /Oq  +  =  c//6  since  />„-  0 

The  remainder  of  the  terms  in  equation  (A  11)  are  unchanged,  thus 

c/V  3  df 
and 

k_  .  f 

i _  # 


k  since  ko  3  1 


where 


f  = 


s; 


zr; 


(A  15) 


from  equation  (A  12)  and  (A  13) 


(A  16) 


substituting  equations  (A  14)  and  (A  15)  into  equation  (A  16) 


<*>x  =  i  -  %Z~  +az*cPx _ cs 

•ci  owe  ~  <ct<» 


Vx 


:<^x 


owe 


<35xXo+  ^ 


or 


<#>x  _ 


<^x  Xa+2k“u 


=  <*x 


agx 


This  is  equation  (17)  on  page  23  .  Since  fuel  depletion  effects  are 
not  considered  during  the  short  simulated  transients,  the  fuel  number 
density  U  is  assumed  constant  and  cX^is  invariant  in  the  time  and  the 
Laplace  domains. 
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APPENDIX  B 

DEVELOPMENT  OF  THE  REACTOR  HEAT  TRANSFER  FUNCTIONS 


The  previously  developed  equations  describing  heat  transfer  from  the 
reactor  to  the  coolant  are 


where 


and 


Tp(*>  *  i±  —fr1  =  TavIO  +  ^  P(t) 


TAv(t}[l*2  %] 


t  12 


dTayfel 


*  TP  (t)+ 2 (tl 


TavCO  =  Ctl  "1 


(B  1) 

(B  2) 
(B  3) 


Tp(t)  ■  Average  core  material  temperature  (°F) 
TAV(t)  ■  Average  reactor  coolant  temperature  (°F) 

Tho ( t )  =  Reactor  coolant  outlet  temperature  (°F) 

Tci (t)  *  Reactor  coolant  Inlet  temperature  (°F) 

P(t)  *  Total  power  generated  In  the  core  (Btu/sec) 


F 

Si 

C 


'FM 


*  Total  thermal  capacity  of  core  materials  (Btu/°F) 

■  Total  thermal  capacity  of  coolant  in  the  core  (Btu/°F) 

*  Specific  heat  of  coolant  (Btu/blm°F) 
a  Coolant  mass  flow  rate  (lbm/sec) 

*  Total  heat  transfer  coefficient  (Btu/sec°F) 


T.  = 

Ti  = 

1z  * 


Ch 


m*C 

C  F 
hpn 

Hpm 


61 


Let 


P(t) 

3  Po  + 

</*P(t) 

TF(t) 

aTf0  + 

<ArF(t) 

AV(f) 

3  TAV  + 

<AAV(t) 

*><*> 

Aio(t> 

(B  4) 


Tc1 (t)  =  Tcio  +  J* Tc^ (t) 

where  the  lowest  zero  subscript  denotes  a  steady  state  value  and 
the  delta  prefix  a  small  perturbation  about  this  value. 

From  substitution  of  equations  (B  4)  into  equations  (B  1),  (B  2), 
and  (B  3),  with  the  reactor  in  a  steady  state  at  ti  0, 


Tfq  a  Tav,,  ^  Po 

Tav.<H-2  =■  7>„  +  2  Tc,.  (B  5) 

Tav„  =  -j-  (1h.,  »  1c,.) 


and 

(B  6) 

a  JTavIo)  -  e/%o( o)  a  c/*7c*  (o)  sr  O 

Impose  the  perturbations  at  t*0.  Substituting  equations  (B  4)  and 
(B  5)  into  equations  (B  1),  (B  2),  and  (B  3)  with  the  initial  conditions 
given  in  equations  (B  6), 


</TpCO  *•  i±  it  53  c/’PCt)  (B  7) 

</lAv(tl  [*  +  z%]  +  1*  a  o^F(t)  +  2  ~  o^Tc^*)  (B  8) 

</TAV(tl  *  I  prHo  it)  t  <fTcK  (t)]  (B  9) 
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Taking  the  Laplace  transform  of  equations  (B  7),  (B  8),  and  (B  9) 
with  the  initial  conditions  given  by  equations  (B  6) 


Tf(3)  [sfl+i]  =  Tav(S)  +  ^  PCS, 

(B  10) 

Tav(s)  [slz  +  i+n  ±|"]  »  Tf(s)  +  z  ^  7c(s) 

(B  11) 

TaVCS)  =  \  [  Tw  (S)  ¥  Tc  <S)] 

(B  12) 

where  the  delta  prefix  and  lowest  subscript  have 

been  deleted  for 

readability 

Substituting  equation  (B  12)  into  equations  (B  10)  and  (B  11) 


(B  13) 
(B  14) 


Tr  (slin)  =  ( T«  *  Tc)  + 

|(TH*Te)(s1i+i  +  »lj|)  =  + 

where  the  s  domain  dependence  notation  has  been  deleted  for 
readability. 

Solving  equation  (B  13)  for  Tp  and  substituting  into  equation  (B  14) 

|(Tm  t  Tc)  (sTittti  %)(^n)  =  \  (Tu+Tc)  (3ii+ 1)  Tc 

After  expanding  the  products  and  collecting  terms 

TH[^(^4tT^)3+|]=-TC[s««+(a4-:f)s-IH,> 

After  multiplying  through  by  and  solving  for  TH(  (s)) 

o 
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and  all  notation  has  been  restored. 


where  %  — 

^  12 

This  is  equation  (26)  on  page  28. 


APPENDIX  C 


DEVELOPMENT  OF  THE  HEAT  EXCHANGER  TRANSFER  FUNCTION 


The  previously  derived  heat  exchanger  heat  transfer  equations  are 
2Ts(0-  TmoOtl-K*! 


Tco<t)  = 


[t  +  Kl] 


£l  +  TSU)  =  T  [  THi(t)  +  lc0(0  ]  +  Pi.lt) 

Wtm 


d* 


where  THi-(t)  3  Heat  exchanger  coolant  inlet  (°F) 

TC()(t)  »  Heat  exchanger  coolant  outlet  (°F) 

T^(t)  =  Saturated  steam  pressure  ( PS I Q ) 

P^(t)  =  Power  delivered  by  heat  exchanger  (Btu/sec) 

Cs  =  Total  thermal  capacity  of  heat  exchanger 

metal  and  secondary  water  and  steam  (Btu/°F) 


(C  1) 

(C  2) 


and 


C  s  Specific  heat  of  coolant  (Btu/Ibm°F) 
hTM  =  Total  heat  transfer  coefficient  (Btu/°F.sec) 
Coolant  mass  flow  rate  (Ibm/sec) 


*1  * 


zmHC 

Htm 


let  rH1(t)  •  rH1 

0 

+ 

TC0<1)  ’  TC 

+ 

«**,<*> 

(C  3) 

Ts(t)  =TS 

+ 

A><‘> 

P,(t)  =  P. 

L 

+ 

JPL(t) 

where  the  lowest  zero  subscript  denotes  a  steady  state  value  and  the 
delta  prefix  a  small  perturbation  about  this  value. 
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With  the  plant  in  an  initial  steady  state  at  t  <  0,  equations  (C  1) 


and  (C  2)  are 


Tc„  - 


2T50  “  [  1+  <1  1 


[l  *  K,} 

Ts0  =  k  (T^.tTe.„) 
<^K0{p)  -  O 

(PTu^(o)  =  cPTq0[o)  =  o^Ts(o)  =  o 


(C  4) 


(C  5) 


(C  6) 


Impose  the  perturbations  substituting  equations  (C  3),  (C  4),  and 
(C  5)  into  equations  (C  1)  and  (C  2) 


2<P Ts-c/^Tmj  (1  -  *i) 

<*Tc°~  (1+  Ki) 


(C  7) 


&  £  fhM  ♦  A(t)  =  i  Ac„w]  - 


(C  3) 


Taking  the  Laplace  transforms  of  equations  (C  7)  and  (C  3)  with  the 
initial  conditions  given  by  equation  (C  6) 


JTsC  s)  ^ 


where  Yc  -  Cs 
Htm 


I  [</Yc0(s)  +  (pT».  (S)  ]  -  y^o^LiS) 


I5S  +  1 


</TCo(9)  =  (2</Ts(s)-o^(s)  [1-lei]} 

These  are  equations  (34)  and  (35)  on  page  33. 
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APPENDIX  D 


DEVELOPMENT  OF  THE  PRIMARY  PIPING  TRANSFER  FUNCTIONS 

1 .  Transport  Delay  Transfer  Functions 

The  transport  delay  of  the  coolant  between  the  outlet  of  the  reactor 
and  the  heat  exchanger  inlet  plenum  can  be  expressed  as 

THi  (t>  '  O  ’> 

P 

where  T„.  (t)  =  Heat  exchanger  inlet  plenum  coolant 
p  temperature  (°F) 

7^  (t)  =  Reactor  coolant  outlet  temperature  (°F) 

Y3  =  Transport  time  delay  (sec) 

Rearranging  terms  and  expanding  equation  (D  1)  in  a  Taylor  series 

Th„<M  =  +  4  ...  (D  2) 

For  slow  temperature  changes,  the  second  order  and  higher  terms 
in  equation  (D  2)  may  be  ignored, 

(D3) 

Let 

THi(t)  =  Tm  +  </th  (t)  (D  4) 

po  ip 

THo<‘>  •  Th  *  <AHo(t> 

00 

where  the  lowest  zero  subscript  denotes  a  steady  state  value  and 
the  delta  prefix  a  small  perturbation  about  this  value. 
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From  substitution  of  equations  (D  4)  into  equation  (D  3)  and  for 


the  temperatures  in  a  steady  state  at  t^O 


and 


=  c/THo(t)  =  o 


(D 

(D 


Impose  the  perturbations  at  t=0.  Substituting  equations  (D  4)  and 
(D  5)  into  equation  (D  3)  with  the  initial  conditions  given  in  equation 
(D  6) 

/TH.p(t)+  ?3  ^  /nj.p£t)  =5  c/'THgCt)  (D 


Taking  the  Laplace  transform  of  equation  (D  7)  and  solving  for 


Ho 


(PTh^ls)  _  1 

<^THo(S)  X  +  73S 


This  is  equation  (39}  on  page  36. 

By  following  a  similar  development,  the  transport  delay  between  a 
perturbation  in  the  heat  exchanger  coolant  outlet  temperature  and  the 
resulting  perturbation  in  the  reactor  inlet  plenum  coolant  temperature 
can  be  expressed  by  the  transfer  function 

cflcoCs)  1  +•  TiS 

where  =  Reactor  inlet  plenum  coolant  temperature 

perturbation  (°F) 

c/'Tco  =  Heat  exchanger  coolant  outlet  temperature 
perturbation  (°F) 

I.  =  Transport  time  delay  (sec) 


This  is  equation  (40)  on  page 


2,  Mixing  Delay  Transfer  Function 

Let  the  reactor  and  heat  exchanger  coolant  inlet  plenums  be  represent¬ 
ed  by  the  control  volume  shown  in  Figure  D-l. 

M  =  Mass  of  coolant  in  plenum  (Ibm/sec) 

T  =  Average  plenum  coolant  temperature 
(°F) 

T.  =  Plenum  inlet  coolant  temperature 
1  (°F) 

T  =  Plenum  outlet  coolant  temperature 
0  CF> 

m  =  Coolant  mass  flow  rate  (Ibm/sec) 

Figure  D-l.  Mixing  Volume 


m,T0 


If  perfect  mixing  is  assumed,  the  coolant  outlet  temperature  is  equal 
to  the  average  coolant  temperature  in  the  plenum.  At  heat  balance  on 
the  volume  requires  that  the  net  heat  flowing  into  the  plenum  be  equal 
to  the  increase  in  stored  energy.  Assuming  no  ambient  heat  losses,  and 
assuming  the  temperature  of  the  plenum  structural  material  remains 
constant,  the  heat  balance  per  unit  time  over  a  time  interval  At  is 

m AtC  [T,(t)-T.(t)]  =  MC  -  mAtC^j^ 


where  C  =  specific  heat  of  coolant 
The  resulting  equation  is 


m 


dfc 


+To{t\  -T\  (t) 
dt 


(D  8) 
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where  T  ®  —  =  mixing  time  delay, 

m 


V‘>  *  T,0  +  /Tf(t) 
T0(t)  =  Too  +  A„(t) 


(D  9) 


where  the  lowest  zero  subscript  denotes  an  initial  steady  state  value 
and  the  delta  prefix  a  small  perturbation  about  this  value. 

Substituting  equations  (0  8)  into  equation  (D  9) 

1  ^  +  JYoCt)  =  (DIO) 


Taking  the  Laplace  transform  of  equation  (DIO)  and  solving  for 

A' 


This  is  equation  (41)  on  page  37, 
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PWRSG  INSTANTANEOUS  POWER  DELIVERED  TO  TURBINE  < 

PWRSGT  INSTANTANEOUS  »DWER  DELIVERED  TO  THE  TURBINE  (BTJ/SEC) 
A  INSTANTANEOUS  TJRBINE  THROTTLE  VALVE  SETTING 
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